Predicting Permissive Mutations That Improve the Fitness of A(H1N1)pdm09 Viruses Bearing the H275Y Neuraminidase Substitution

ABSTRACT Oseltamivir-resistant influenza viruses arise due to amino acid mutations in key residues of the viral neuraminidase (NA). These changes often come at a fitness cost; however, it is known that permissive mutations in the viral NA can overcome this cost. This result was observed in former seasonal A(H1N1) viruses in 2007 which expressed the H275Y substitution (N1 numbering) with no apparent fitness cost and lead to widespread oseltamivir resistance. Therefore, this study aims to predict permissive mutations that may similarly enable fit H275Y variants to arise in currently circulating A(H1N1)pdm09 viruses. The first approach in this study utilized in silico analyses to predict potentially permissive mutations. The second approach involved the generation of a virus library which encompassed all possible NA mutations while keeping H275Y fixed. Fit variants were then selected by serially passaging the virus library either through ferrets by transmission or passaging once in vitro. The fitness impact of selected substitutions was further evaluated experimentally. The computational approach predicted three candidate permissive NA mutations which, in combination with each other, restored the replicative fitness of an H275Y variant. The second approach identified a stringent bottleneck during transmission between ferrets; however, three further substitutions were identified which may improve transmissibility. A comparison of fit H275Y variants in vitro and in experimentally infected animals showed a statistically significant correlation in the variants that were positively selected. Overall, this study provides valuable tools and insights into potential permissive mutations that may facilitate the emergence of a fit H275Y A(H1N1)pdm09 variant. IMPORTANCE Oseltamivir (Tamiflu) is the most widely used antiviral for the treatment of influenza infections. Therefore, resistance to oseltamivir is a public health concern. This study is important as it explores the different evolutionary pathways available to current circulating influenza viruses that may lead to widespread oseltamivir resistance. Specifically, this study develops valuable experimental and computational tools to evaluate the fitness landscape of circulating A(H1N1)pmd09 influenza viruses bearing the H275Y mutation. The H275Y substitution is most commonly reported to confer oseltamivir resistance but also leads to loss of virus replication and transmission fitness, which limits its spread. However, it is known from previous influenza seasons that influenza viruses can evolve to overcome this loss of fitness. Therefore, this study aims to prospectively predict how contemporary A(H1N1)pmd09 influenza viruses may evolve to overcome the fitness cost of bearing the H275Y NA substitution, which could result in widespread oseltamivir resistance.

The in silico effect of each substitution on protein stability was calculated by FoldX, where a negative number represents increased protein stability and a positive number represents reduced stability. For comparison, the effect of H275Y in reducing protein stability was calculated and was shown to reduce in silico protein stability more substantially than any other substitution. The candidate substitutions S286G, S299A and S95N that were selected for experimental analysis are indicated by blue dots (Panel A) or highlighted in blue (Panel B). The permissive substitutions previously evaluated by Butler et al. [1] are indicated by red dots (Panel A) or highlighted in red (Panel B). Of note the V106I substitutions was not chosen for further evaluation as it has been replaced with the S200N substitution in currently circulating strains.  The PCR cycling conditions were as follows: N1-Rev at 4.5µM and 3µL water. The PCR cycling conditions were the same as described above, but using 20 cycles instead of 7. The joined products were purified, diluted to 3 ng/µL and used as a template for a second round of mutagenesis and joining PCR reactions. The joined products from the second round were then used for ligation.

Gibson assembly for preparing plasmid library
The NA PCR library was utilised to generate an NA plasmid library via Gibson assembly and high efficiency transformation. The pHW2000 plasmid was linearized by BsmBI enzyme digestion and dephosphorylated with Thermosensitive Alkaline Phosphatase (Promega, USA).
The assembly for the joined products were set up in duplicates, where 100 ng of purified vector were also plated out for colony counting. At least five independent transformations were carried out for each preparation of the assembled products. The following day there were between 200,000 to 400,000 colonies on each agar plate, as determined by colony counting from the diluted agar plates. After pooling colonies from the undiluted agar plates (2,000,000 x 5 = 10 6 colonies total), they were cultured for 4 hours at 37ºC and then maxi-preps were performed using the EndoFree® Plasmid Maxi kit (Qiagen, Germany).
To gain insight regarding the distribution of codon mutations in the plasmid libraries, a total of 54 clones were picked from the diluted agar plates (with equivalent numbers of clones taken from the three replicates), and Sanger sequencing was performed. The distribution of codon mutations in the Sanger-sequenced clones were further analysed using a custom python script (https://github.com/jbloomlab/SangerMutantLibraryAnalysis).The three plasmid libraries (i, ii, and ii) were deep sequenced to determine if all possible substitutions had been comprehensively sampled.

Reverse genetics for generation of virus library
Each NA plasmid library was used for reverse genetics to generate a virus library ( Figure 1C) [2,4,5]. Four ferrets were experimentally inoculated with 500 µL containing 10 4.7 TCID50 of pooled virus library (day 0), as previously described [10]. One ferret was experimentally inoculated with the SA16-H275Y virus as a control. Each, experimentally infected ferret was then co-housed with a naïve contact recipient (direct contact 1) 24 hours post-inoculation. Nasal washes were performed daily on direct contact 1 ferrets and nasal wash samples were analysed for infection by qPCR [11]. On the first day that nasal wash samples from direct contact 1 ferrets were qPCR positive for influenza virus, the animal was removed from the cage, and co-housed with a second naïve recipient (direct contact 2). Similarly, nasal wash samples from direct contact 2 were monitored for influenza virus. On the first day that nasal wash samples from direct contact 2 ferrets were qPCR positive for influenza virus, these animals were placed in aerosol cages, adjacent to a third set of naïve recipients (aerosol contacts). Due to limited animal numbers, the SA16-H275Y virus was only passaged once through ferrets (Experimentally infected animals to Direct Contact 1).
In the experiments described, ferrets were nasal washed every day and weight and body temperatures were collected as previously described [12]. Experimentally infected ferrets were euthanized on day 4 of the experiment, and all other animals were euthanized on day 14 of the experiment. Viral titres in nasal wash samples were determined by qPCR [11] and TCID50 assay [8]. for sequencing. The full genome sequencing was done after amplification of all genes using primers previously described [14]. Sequencing of amplified PCR products were done at the Australian Genome Research Facility, on the HiSeq 2500 platform (2x 150 PE reads, 15 million reads per sample).

Analysis of deep sequencing data and bioinformatics
The and ratio of synonymous to non-synonymous mutations in ferret nasal wash samples was calculating by measuring π and πS/πN using the SNPgenie software [17]. The nucleotide mutation frequencies in donor:recipient pairs from the eight contact transmission pairs and for aerosol transmission pairs were also used to estimate transmission bottleneck sizes using the beta-binomial sampling method developed by Leonard et al [18]. This statistical method takes the stochastic dynamics of viral replication in recipients into account and further considers variant calling thresholds. For our analysis, a minimum variant calling threshold of 1% was utilised to estimate bottleneck size to include a greater number of sites, as was done by Poon et al. in a human household transmission study [19]. A more conservative estimate of the bottleneck size was also calculated, using a minimum variant calling threshold of 3% similar to Leonard et al. [18].
The dms_tools2 and mapmuts pipelines developed within the laboratory of Dr. Jesse Bloom (Fred Hutchinson Cancer Research Center, Washington, USA) were also used to produce codon-count files, which were then used for amino acid preference analysis [1-3, 6, 15]. Amino acid preferences were used to calculate Shannon entropy, which is a measure of mutational tolerance of the NA protein in experimental settings. Shannon entropy was also calculated for all natural sequences with the H275Y mutation in the GISAID database, after sequence clean-up and alignment using MAFFT [20]. Pymol (Schrödinger, USA) [21]